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Abstract We study the parametric wave interaction in qua- 
dratic nonlinear media with randomized distribution of the 
ferroelectric domains. In particular, we discuss properties of 
second and cascaded third harmonic generation. We derive 
analytical formulas describing emission properties of the sec- 
ond and third harmonics in the presence of domain disorder 
and show that the latter process is governed by the character- 
istics of the constituent processes, i.e. second harmonic gen- 
eration and sum frequency mixing. We demonstrate the role 
of randomness on various second and third harmonic gener- 
ation regimes such as Raman-Nath and Cerenkov nonlinear 
diffraction. We show that the randomness-induced incoher- 
ence in the wave interaction leads to deterioration of conver- 
sion efficiency and angular spreading of harmonic generated 
in the processes relying on transverse phase matching such 
as Raman-Nath. On the other hand forward and Cerenkov 
frequency generation are basically insensitive to the domain 
randomness. 



1 Introduction 

Periodic poling has been commonly used to realize efficient 
frequency conversion in ferroelectric crystals with quadratic 
nonlinearity d . The spatially periodic modulation of the sign 
of x^^^ nonlinearity ensures the synchronization of the phases 
of interacting waves via the so called quasi phase match- 
ing EKSllll. In the simplest case of ID periodically polled 
structure [Fig.l(a-b)] and collinear second harmonic gener- 
ation the period of nonlinearity modulation A satisfies the 
following relation |q| = 27r/A = k2 — 2ki. Here q is the 
reciprocal vector of the nolinearity modulation and ki , /c2 are 
wave vectors of the fundamental and second harmonic, re- 
spectively. The quasi-phase matching technique is so versa- 
tile that by modulating spatially the domain pattern one can 
also spatially shape the wavefront of generated waves as has 



been demonstrated in case of Bessel and Airy beams SHG ||5j 
01, however, because of resonant character of the quasi phase 
matching, perfectly periodic structure is efficient only for a 
particular choice of the wavelengths of interacting waves. 
Few methods have been proposed in order to extend the oper- 
ational bandwidth of the QPM technique. They all rely on en- 
gineered poled structures and involve, for instance, structures 
consisting of periodically poled sections with different pe- 
riods or multi-period and quasi-periodic structures, i.e. non- 
linear superlattices ir7l[8ll9l [TQl[TT1l . Another possibility offer 




Fig. 1 (color online) (a) Illustrating the concept of periodically 
poled quadratic crystal. Arrow indicates the reciprocal grating vec- 
tor, (b) Domain structure in commercial sample of periodically 
poled lithium niobate crystal with nominal poling period of 14/im 
(visualized via the second harmonic nonlinear microscopy), (c-d) 
Schematic representation of the quadratic crystal with (c) one- and 
(d) two-dimensional random spatial domain distribution. Arrows 
represent various reciprocal grating vectors, (e) Random ferroelec- 
tric domain pattern in SBN crystal (visualized via SH nonlinear mi- 
croscopy (V2\ ). 'C in (b) and (e) indicates the direction of the opti- 
cal C-axis. 
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Fig. 2 (color online) Illustrating the phase matching conditions for 
the second and third harmonic generation via nonlinear diffraction in 
quadratic medium with randomized ID ferroelectric domain struc- 
ture. /c2 and ks represent wave vectors of the second and third har- 
monics, respectively. Green and blue rings (with corresponding radii 
k2 and /cs) define all possible emission directions of the second and 
third harmonics. 



media with random or randomized distribution of the sign of 
nonlinearity. In fact, Horowitz et al f\3l and Kawai et al |[T4l 
were first to demonstrate broadband parametric interaction in 
unpoled strontium barium niobate (SBN) crystals. As-grown 
SBN crystals [Fig.l(d-e)] consist of elongated ferroelectric 
domains oriented along its optical (C) axis having random 
distribution of their size and orientation flSl . As such they 
represent disordered nonlinear medium. Such media have at- 
tracted a considerable attention in recent years due to their 
potential benefits for nonlinear parametric processes II16II17I 
[T8]| . In the particular case of SBN crystal its spatially ran- 
dom distribution of nonlinearity [see Fig. 1(d)] provides phase 
matching for second-order nonlinear processes for practically 
any incident wavelength [191 in its transmission band. The 
extensive studies of the second-order nonlinear processes in 
random SBN crystals resulted in demonstration of a num- 
ber of effects including e.g., angular broadening of the emit- 
ted harmonics II20II21L multi-frequency conversion II22II23I 
l24l . conical and planar frequency emission II25[|26[|27]| and 
its application for optical autocorrelator II28II29I . Few recent 
works have demonstrated beneficial effect of controlled ran- 
domization of originally periodic domain pattern on second 
harmonic generation II30II31I It should be also mentioned that 
the randomness of the ferroelectric domains pattern can be 
introduced inadvertently in the fabrication of QPM structures 
via periodic poling [see Fig. 1(b)]. This typically has an ad- 
verse consequences on the intended performance of the QPM 
device leading, e.g. to decreased frequency conversion effi- 

It has been also shown recently that under certain condi- 
tions the ferroelectric crystal with random domain distribu- 
tion can be used for broadband cascaded third harmonic gen- 
eration |[T9ll36[|37[|38[[39l . In such case the randomness of the 
nonlinearity contributes to both constituent processes: second 
harmonic and sum frequency generation. Subsequently, as the 
experiments showed, the spatial light intensity distribution of 
the emitted third and second harmonics are significantly dif- 
ferent. So far, no theoretical analysis of the spatial proper- 



ties of light emissions in the cascaded process has been con- 
ducted. 

In this work we present extensive theoretical studies of 
the non-collinear second and third harmonic emission in fer- 
roelectric crystal with randomized domain structure. Using 
the exact analytical formulas for the intensity of generated 
waves obtained in the regime of strong fundamental beam we 
compare the emission characteristics of both harmonics as a 
function of the degree of disorder and the wavelength. 



2 Sum Frequency Mixing 

Let us consider interaction of two optical waves with ampli- 
tudes Ea , with the frequencies uOa and cj^ in quadratic 
crystal with the spatially modulated sign of nonlinear coef- 
ficient M(r). The nonlinearity will lead, among others, to 
the generation of the third wave E^p at the frequency cJsf = 
Ua + whose amplitude can be expressed as |40|: 



E,M)(xidEaEl j j^M{v)e'^^dv (1) 

where d is the nonlinearity coefficient and integration is per- 
formed over the whole area A of the nonlinear crystal. Here 
q = ka — with ki = ni27r/Ai representing wave vector 
and rii is refractive index of the corresponding wave. In this 
paper we consider periodic structure with modulation only in 
X direction (Fig. [2]) such that: 



N-l 
3=0 



(2) 



where Xj is a position of j-th domain wall, N is number of 
domains and nx.^x.^^{x) = 1 for Xj < x < Xj+i and 
otherwise. In this case integral ([B can be separated into x 
and y directions. Additionally, because M(r) = outside 
crystal integration can be extended into infinity and one gets: 

£;sp(q) oc idEaE; f ' e^^^^d^ / ^ M{x)e''^-''dx 
Jo Jo 

Qy^y 



(X idEaE^LySinc 



where L^^ Ly are the lengths of the crystal in x, y directions 
and M{q^) denotes a Fourier transform of spatially modu- 
lated nonlinear coefficient: 

N 

(4) 

In general the last two terms in the bracket can be neglected (H 
what results in: 



. N 



j=l k=l 



(5) 
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where Ik = Xk — Xk-i is the width of k-\h domain. If the 
domains' widths Ik are not equal but are subject to some ran- 
dom process then one has to take into account the statisti- 
cal properties of this process. Additionally when generation 
results from many domains the intensity of generated wave 
/sF = |£^sFp is represented by an average over different real- 
ization (< • > denotes averaging over domains width /): 



(/sp) oc d^Lpahsinc 



<|M(gJp). (6) 



The expression (|M(g)p) has been discussed before Rll4T1l. 
We adopted here the approach as described by Le Grand et ah 
with some modifications that allow us to extend the results of 
SHG to SFM and subsequently to the cascaded THG. 

If domains width varies randomly it can be described by 
an arbitrary Probability Density Function (PDF) w(l) with 
characteristic function (^{q) = J w{l) ex.p{iql)dl and for \(p{q)\ 
1 and large number of domains (N 1) one obtains: 



(|M(g)|2) = ^(^l-2i?e 



f{q) + 1 



AN 



(7) 

what is exactly the formula derived by Le Grand et al. BTIl . 

The effect of averaging and above discussed simplifica- 
tion on SH intensity is shown in the Fig. [3] This picture com- 
pares the square modulus of Fourier coefficients calculated 
for three different cases: the gray solid line depicts |M(g)p 
with M{q) calculated from Eq. Q, red dashed line shows 
value of (|M(g)p) where \M{q)\'^ is calculated as for grey 
line and then averaged over 300 randomly generated sam- 
ples, finally the dotted blue line shows (|M(g')p) calculated 
according to simplified Eq. d?]). One can see that results ob- 
tained with the simplified equation are consisted with numer- 
ically averaged SH intensity calculated using the full model 
for both very small (a) (cr=0.05/im) and large (b) dispersion 
(a=lOfim). All curves ware calculated assuming number of 
domains to be N=500. 

In order to calculate the averaged SF intensity (/sf) one 
has to choose a function that describes statistical distribution 
of the domain widths. A natural choice would be Gaussian 
distribution. However Gaussian model is only applicable to 
rather small values of dispersions a (in comparison to the av- 
erage domain size Im) because for dispersions comparable or 
bigger then Im it results in negative values of domain width. 
The more appropriate choice is therefore the Gamma distri- 
bution which is defined for / G (0, oo) with the PDF as 



w{l) 



k-lL 



(8) 



where k and 6 are the so called shape and scale parameters, 
respectively and F is gamma function. In this case the charac- 
teristic function is (^{q) = (1 — iqO)~^ . Maximum of the dis- 
tribution, i.e. the most likely value (or mode) is = {k — l)0 
while standard deviation is a = Vk6^. Gamma, unlike the 
Gaussian distribution, enables one to account for all possible 
dispersions regardless of the mode value. 
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Fig. 3 (color online) Fourier coefficient of the sum frequency mix- 
ing process as a function of the phase mismatch parameter q for a 
periodic cr=0.05/im (a) and randomized cr= lO/xm (b) domain pat- 
tern. Solid (grey) and dashed (red) lines represent formula Eq.(|4]) 
^ calculated for particular realization of the domain distribution and 
averaged over 300 numerically generated random structures, respec- 
tively. Dotted (blue) line represents simplified formula valid for 
large number of domain. The blue line almost exactly overlaps with 
the dashed red line. 



In this paper, for the sake of clarity of the presented re- 
sults, in all calculations we use Gamma distribution with fixed 
mode value Im = IjJ^m and we will characterize the random- 
ness of the structure by its standard deviation a. 



3 Second Harmonic Generation 

We consider the situation where the fundamental beam with 
the wavelength A is incident onto the sample along the y di- 
rection (Fig. O. Two photons of the fundamental harmonic 
2/ci, generate non-collinear SH with wave vector /c2 emitted 
at the angle a2 such that the corresponding phase mismatch 
equals to q2 = k2 — 2ki (triangle OAB in Fig.O. According 
to the formula Eq. ([3]) (with Ea = Ei^ = Ei) the amplitude 
of the generated SH electric E2 field is given as 



£^2(^2) oc idLy Iisinc 



Q2yLy 



(9) 



and the averaged SH intensity {I2)'. 



{h) oc d'ltKsinc 



(|M(g2j|'). (10) 



This formula relates the stochastic properties of the fer- 
roelectric domain distribution to the spatial average intensity 
distribution of the SH. The most universal way to analyze this 
result is to use a 2D map plotting the Eq. (fTOb as a function 
of arbitrary phase mismatch components q^ and qy. This is 
demonstrated in Fig. |4] where we show the average SH in- 
tensity ((/2)) calculated for a sample with nominal domain 
width /^=l/im (what corresponds to poling period A=2jj.m) 
and four different values od dispersion a: plots (a-d) corre- 
spond to the increasing degree of the randomness ranging 
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Fig. 4 (color online) Illustrating the effect of randomness on the 
strength of generation of the second harmonic in periodically poled 
media with random perturbation of the regular domain structure. 
Shown is the 2D map of the averaged intensity (in logarithmic scale) 
of the SH (EqlTot as a function of the phase mismatch and Qy. 
In all graphs the average domain size is I = 1/im while the dis- 
persion of the domain size varies from cr=0.05/i m to cr= lO/xm. The 
dashed and dotted white curves represent all possible directions of 
the wave vectors of the second harmonics of the fundamental waves 
of A = 1.2/im and A = 1.73/im, respectively. The inset depicts 
details of the map in the vicinity of (qx = 0,qy = 0). Periodic os- 
cillation of the intensity reflect the Makers fringes of the forward 



from cr=0.05/im (almost ideal periodic structure), cr=0.5/im, 
cr=l/im and cr=10/im (strongly disordered structure). In case 
of periodic structure Fig. |4] (a) there are five characteristic 
lines/traces which correspond to the strong Fourier compo- 
nents of M{qx^qy). These are: i) a vertical line at = 
which represents the Cerenkov SH radiation I42II43I (marked 
as Ce), ii) two symmetrically positioned horizontal lines at 
g^=zb3.1/im~^ corresponding to first order Raman-Nath res- 
onance I44II45I (marked as RN±i), iii) two symmetrically po- 
sitioned horizontal lines at q^=±9Ajj.m~^ corresponding to 
the third order Raman-Nath resonances (marked as RN^s). 
When the dispersion of the domain size increases [(Fig.|4](b)] 
the Raman-Nath emission lines shift towards the center and, 
at the same time, become broader and weaker until they dis- 
appear for strongly disordered structure (Fig.|4l(c-d)). On the 
other hand the Cerenkov SH signal can either increase or de- 
crease depending on the value of q^ . 

The dashed and dotted white curves represent phase mis- 
match (q2x^Q2y) calculated for two arbitrary chosen wave- 
lengths of the fundamental wave: X=l.2fi (dashed line), and 
A=1.73/im (dotted line) using reported refractive index of 
LiNbOs crystal 1461 . If we fix fundamental beam wavelength 
(say, A=l .2/im) it is clear that for almost periodic domain dis- 
tribution the SH is efficiently generated only in those (qx^ Qy) 



Fig. 5 (color online) Normalized intensity of the SH signal as a 
function of emission angle a2 and wavelength of the fundamental 
beam A for various degree of disorder: (a) cr=0.05/im, (b) cr=0.2/im, 
(c) cr= 1/im and (d) cr= lO/xm. 



points where the phase mismatch curve intersects the lines 
corresponding to large value of the Fourier spectrum of the 
nonlinearity modulation. Those characteristic points result in 
a wavelength dependent spatial emission pattern of SH. 

Fig. [5] shows the intensity of the SH (Eq. ([TOb ) as a func- 
tion of fundamental wavelength and SH emission angle 0^2. 
For a periodic structure [Fig.[3a)] SH intensity exhibits clear 
peaks where the emission angles 0^2 are defined by partial or 
full phase matching conditions. For Cerenkov SH emission, 
where ^^2^ = 0, the emission angle can be calculated as: 



cos(a^^) 



ni 
^2' 



(11) 



and for m-th order Raman-Nath where q2x = m27T/A, emis- 
sion angle is: 

A 



sin(ar-) 



2An2 



(12) 



Those wavelength dependent emission angles are plotted in 
the Fig|5] (a) as a red solid line (marked as Ce) for and 
as red dashed lines (marked as RNi, RN3) for af^"". As it 
has been mentioned earlier, it is clearly seen that the emit- 
ted Raman-Nath SH exhibits enhanced angular spreading for 
strong domain dispersion while Cerenkov emissions is less 
sensitive to disorder. 

The horizontal dashed and dotted lines in Fig. [5] point to 
two arbitrary chosen wavelengths of the fundamental beam 
A=1.2/im (white dashed line) and A=1.73/im (white dotted 
line). The corresponding second harmonic angular intensity 
profiles for cr=0.05/im (red solid lines) and cr=10/im (blue 
dashed line) are depicted in FiglQ The strong modification 
of the emission pattern with high degree of randomness is 
evident, however it affects SH in different ways depending 
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Fig. 6 (color online) SH intensity as a function of the emission an- 
gle a2 for cr=0.05/im (cross section along dotted line in Fig. [3 a)) 
and cr=10/im (cross-section along dotted line in Fig.[5td)). Funda- 
mental wavelength A=L73/im (a) and A=1.2/im (b). 



subject to intensity modulation. The most efficient SH emis- 
sion is expected where the Raman-Nath and Cerenkov signals 
overlap (as it happens for X=l.l3fim and A=0.96/im) what 
constitutes the regime of so called nonlinear Bragg diffrac- 
tion |47l, i.e. SH emission under simultaneous fulfillment of 
the transverse and longitudinal phase matching. For strong 
disorder the peak intensity deteriorates due to decrease in RN 
intensity. However at some wavelengths (as for X=l.2fim) the 
Ce intensity increases. This strong dependence of the Cerenkov 
signal on the wavelength is a result of the interference of 
many SH components emitted by each ferroelectric domain BSl . 
When the initially periodic domain pattern becomes disor- 
dered the interference effect gets weaker and the total Cerenkov 
signal is determined by the sum of the intensities of contribut- 
ing waves. In effect, the Cerenkov signal which was strong 
for ceratin fundamental wavelength (due to a constructive in- 
terference) is weakened while the one which was initially 
weak (due to destructive interference) increases. 



on the wavelength and the type of interaction. For A=1.73/im 
[FiglSa)] we can see that increasing of a suppresses intensity 
of the Cerenkov (Ce) signal. On the other hand for A=1.2/im 
[Fig 13b)] Ce emission significantly increases while first and 
third order Raman-Nath peaks (RNi and RN3) practically 
vanish. For both wavelengths strong dispersion produces zero 
order Raman-Nath (RNq) reciprocal vectors. 
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Fig. 7 (color online) Intensity of the Cerenkov and Raman-Nath 
second harmonic signals as a function of fundamental wavelength 
for different dispersions of the domain size: (a) cr=0.05/im and (b) 
cr=10/im. Red solid line (Ce) represents Cerenkov intensity, blue 
dashed line (RNi) is the first order Raman-Nath and the green 
dashed-dotted line (RN3) represents the third order Raman-Nath. 



Fig. [7] shows Cerenkov SH (red solid line), Raman-Nath 
first (blue dashed line) and third order (green dashed-dotted 
line) SH as a function of wavelength for two different disper- 
sions (j=0.05/im [FigOa)] and a=lOfim [Fig. Ob)]. Those 
plots can be also understood as profiles along vertical red 
curves marked on a the Fig. Oa) as Ce, RNi and RN3 re- 
spectively. One can notice that Cerenkov radiation (Ce) is a 
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Fig. 8 (color online) Intensity of the Cerenkov (Ce), zero (RNo) and 
first (RNi) order Raman-Nath signals as a function of the domain 
dispersion a, for A = 1.2/im and A = 1.73/im (as depicted in 
brackets). 



This is even clearer in the Fig. [8] where we plot the Cerenkov, 
zero and first order Raman-Nath as a function of domain dis- 
persion a. We can see that with growing disorder Raman- 
Nath RNi(A=1.2/im) and Cerenkov Ce(A=1.73/im) signals 
decrease and tends to a constant value as the constructive in- 
terference effects associated with periodicity of the domain 
distribution are no longer relevant for large a. On the other 
hand signal Ce(A=1.2/im) increases due to weakening of de- 
structive interference effect. Zero order Raman-Nath RNo(A=l .73/im), 
initially not present for perfect structure ie. a 0, continu- 
ously increases with a as the high degree of disorder pro- 
vides small but nonzero fourier coefficients for q2x=^ which, 
in case of SH, contribute towards the forward emission. 



4 Third Harmonic Generation 

It has been demonstrated recently that nonlinear crystals with 
random domain distribution enable third harmonic genera- 
tion via cascading of the frequency doubling and sum fre- 
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quency mixing |[T9l[36l[38l[39l . In this section we provide the 
first analytical study of such process by extending the ap- 
proach discussed in the previous two sections. We will con- 
sider the same setting of the one-dimensional sample of the 
quadratic crystal with randomized domain distribution being 
illuminated by the strong fundamental beam. The nonlinear 
interaction leads to simultaneous emission of the second and 
third harmonics. To calculate the amplitude of the third har- 
monic we use Eq.([3l) where Ea = E2 is taken from Eq.® 
and El) = Ei. In this case the amplitude of the electric field 
of the third harmonic Es is 



£^3(q2,q3) oc -d' L;EI hsinc ^ 
xM{q2,)M{qs,) 



Q2yLy 



smc 



QSyLy 



(13) 



where q2 and qa represent the phase mismatch of the first and 
the second process of sum frequency mixing, respectively. 
The phase matching diagram of the TH generation via cas- 
cading process is illustrated in Fig. [21 This graph represents 
relevant processes in this multi-wave interaction. First a pho- 
ton of SH is generated in the frequency doubling process with 
the phase mismatch wave vector of q2 and emitted at the an- 
gle a2 as depicted by the triangle OAB. Then this photon 
together with photon of the fundamental wave (BC) forms 
a photon (OF) with the triple frequency of the fundamental 
where TH and phase mismatch vector qs = ka — k2 — ki . 
The TH wave is emitted in the direction determined by the 
angle 0^3 (angle between OF and OA). In order to revel key 
features of this cascaded process we will first calculate the 



average intensity of the TH wave Is 
in such a single process: 



=< \Es\ > generated 



4 r4 7-3 • 2 1 ^'^y-^y \ - 2 

' ' "^"^ ' - ' smc 



{Is) ^U^/f sine ^ 

x<\M{q2,)\'\M{qs,)f> 



(14) 



where: 



<\M{q2,)f\M{qsjf>= 
X {NHiq^MQsJ ^ NO, 
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(15) 



Oo 



where Oi, O2 are function of phase mismatch ^23^, QSx 
characteristic function (p{q) and where M{qx) was taken from 
Eq. ([5]). When 1 and \(p{q)\ < 1 terms containing Oi 

and O2 can be neglected and above simplifies to: 



^4 r4 r3 7\r2.;„.2 / ^'^y^v \ ^^^^2 



{Is) oc Wd^'Ll^N'smc^ 



QSyLy 



yQ2x) [QSx) h 



(16) 



The above equation shows that the resulting intensity of the 
third harmonic is a product of intensities of two constituent 
sum frequency mixing processes where both processes utilize 
the same Fourier space. 
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Fig. 9 (color online) Normalized intensity of the TH signal as a 
function of the wavelength of the fundamental beam and emission 
angle for very weak (cr=0.05/im - left column) and strong (cr= lO/xm 
- right column) domain dispersion. Pairs of graphs (a,b), (c-d) and 
(e-f) correspond to three distinct interaction mechanisms represent- 
ing different type of the SH emission according to Eq. ([161 They 
include the cases of second harmonic emitted as a Cerenkov (a,b), 
as well as first (c-d) and third (e-f) order Raman-Nath waves. 



Since TH is generated in a cascaded process and its strength 
depends on the strength of the SH generated in a first step 
it is natural to analyze how SH generated in different pro- 
cesses like Cerenkov or Raman-Nath affects spatial TH emis- 
sion. Fig.[9ta) show TH intensity {Is) calculated according 
to Eq. (O as a function of the fundamental wavelength and 
TH emission angle 0^3 resulting only from SH generated in 
a Cerenkov-like process, i.e. fixed values of SH phase mis- 
match vector to q2y =0 and q2x = Q2 sin(a2^^). One can 
clearly observe a number of curves of enhanced TH inten- 
sity. Periodic structure exhibits strong resonances in fourier 
space [see FigllJa)] and since both processes, SHG and fol- 
lowing SFM utilize the same Fourier space it is expected to 
obtain high TH intensities for processes resulting from either 
Cerenkov or Raman-Nath resonances [45 ). The emission an- 
gle 0^3 of TH generated in double Cerenkov process (ie when 



QSy =0) can be calculated as: 



cos(a^^^^^) = ^ 
^ ^ ^ ns 



(17) 



Angle 0^3^"^^^ as a function of wavelength A is plotted in Fig. [3a) 
as a dashed red line (marked as Ce) and it exactly overlaps 
with one of the TH signals. Similarly the remaining traces can 
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be identified as interaction of Second Harmonic Cerenkov 
with folowing Sum Frequency Mixing of p-th order Raman- 
Nath (ie. qs^ = p27r/A) where the emission angle 0^3 is: 



. / Ce+RN„ N 

sm(a3 ^) 



Xp 



3An^ 



(18) 



For the sake of clarity we plotted only few traces: RN_i, 
RN_3 and RN_5. The strong horizontal line at A=1.73/im 
and the weaker one at A=0.96/im result from the strong SH 
generation at this wavelengths due to Cerenkov resonance 
with RNi and RN3 as shown in Fig.[7j 

When the domain structure becomes strongly disordered 
[FigHb)] all the previously observed Raman-Nath peaks dis- 
appear. However, as in the case of SH also now the Cerenkov 
emission (red dashed line) remains relatively strong and ad- 
ditionally the zero order Raman-Nath (RNq) is present (white 
dotted line). 

Analogously Figs. [3c,d) show average intensity of the 
third harmonic ({I3)) but this time SH is generated as a first 
order Raman-Nath (m = 1) and Figs.[9je-f) third order (m = 
3) Raman-Nath. The formulas for emission angles in those 
processes are: 



cos(af-^^^) 



3ns 



for TH generated as a Cerenkov (red dashed line) or 

A(m + p) 



cos(a3 ) 



3Ans 



(19) 



(20) 



when TH is generated as the p-th order Raman-Nath (indi- 
cated by white dotted line). 



7=10iam 



10"^' - 



10"^' - 



10" 







Eq.(23) a) 






Eq.(24) 
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Fig. 11 (color online) Normalized intensity of the third harmonic as 
a function of the fundamental wavelength and emission angle in the 
regime of almost periodic (cr=0.05/im) and disordered (cr=10/im) 
ferroelectric domain distribution. Plots (a-d) in the left column are 
obtained assuming incoherent contribution to the total TH intensity 
according to Eq. [23] Plots shown in the right column are obtained 
by assuming coherent model Eq.[22l and by averaging the calculated 
TH intensity pattern over 45 different realizations of the domain dis- 
tribution. 



Fig. 10 (color online) Comparison of THG via coherent process ac- 
cording to Eq.(l2n) (gray lines) and incoherent Eq.(|23ll (red lines) 
for a) small dispersion cr=0.05/im and b) large dispersion cr=10/i m. 
The solid blue line represents coherent contributions to the overall 
THG averaged over an ensemble of 120 randomly generated realiza- 
tions of domain pattern [Eq.|22l. Incident wavelength A=L73/im. 



In experimental reality there is no possibility to limit THG 
to only selected processes and so in general the total inten- 
sity of the third harmonic Is would be a sum of TH gen- 
erated from all possible emitted second harmonics. As de- 
scribed previously TH can be generated in a process OBCF as 
schematically shown in Fig El However, the TH with exactly 
the same wave vector can be also generated in an another pro- 
cess as for example the one indicated by ODEF in the same 
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figure. In order to calculate total TH emitted in a particular 
direction as an integral over all possible SH emission angles 
0^2 has to be calculated: 



rcoh 



J £^3(q2,q3) 



da- 



(21) 



and in case of presence of randomness its averaged version is 
given as 



{ll^')^d^Lp!{\ /^^sinc(^ 



X M{q2^)M{qs^)da, 



' q2yLy\ . / q?>yLy 
' smc ' 



(22) 



where averaging again is performed over all realizations of 
domains widths. Both integrals Eq. ([2T]) and Eq. ([22b cannot 
be simplified and can be evaluated only numerically. How- 
ever our previous experiments in Third Harmonic Genera- 
tion in nonlinear crystals with 3 -dimensional fully random 
domain structure indicate that the emission process is in fact 
incoherent, i.e. the intensity of the TH is given as a superpo- 
sition of intensities generated by all constituent second har- 
monic waves ® ,38|. Therefore for large dispersion of do- 
main sizes it is justified to calculate TH signal assuming in- 
coherent emission 



{IT) 



/TT nTT 
(Is) da2 cx / (/sF(q2)) (^srCqs)) da2 
-TT J —TV 



(23) 

where we use (Is) from Eq. ([T6b . 

In Fig{TOl we show angular intensity profile of the third 
harmonic in case of two domain dispersions: a=0.05fim (a) 
and a=lOfim (b) for fundamental wavelength A=l .73/im. The 
grey noisy line represents the TH intensity based on coherent 
superposition of contributing second harmonic waves, calcu- 
lated for a single particular realization of domain distribution 
(Eq. (EB with TH electric field Es taken from Eq.O). The 
solid blue line depicts coherently calculated intensity (same 
as grey line) but averaged over 120 realizations of randomly 
generated structures [Eq. (l22l)1. Finally, the dashed red line 
shows the intensity profile obtained by assuming an incoher- 
ent interaction according to Eq.(l23l). Both models clearly co- 
incide confirming the incoherent character of third harmonic 
generation in that case. The labels next to the emission peaks 
in graph Fig{TOl identify processes responsible for particu- 
lar emission angle. The labeling convention used here and in 
Fig.[TT]is such that the first symbol refers to the second har- 
monic process while the second to the sum frequency mixing. 
Notice how the high degree of randomness leads to increase 
of the forward emission of TH. This is a drect consequence 
of the disorder-enhanced forward emission of the second har- 
monic as seen in Fig. 5(d) and Fig. 6. 

Fundamental wave A=l .73/im is special case when Raman- 
Nath peaks overlap with Cerenkov. In order to ensure that the 
process is incoherent for broader spectrum we calculated the 
angular TH intensity distribution for wavelength range 0.9- 
2.1/im. The results are shown in Fig. [TT] which depicts the 



angular TH intensity distribution for few values of the disper- 
sions of the domain width, a = 0.05/im, a=0.2fim, a=lfim 
and cr=10/im. Maps in the left column (a-d) are obtained as- 
suming incoherent contributions to the total intensity. On the 
other hand, maps (e-h) represent coherent superposition of 
constituent harmonic waves (averaged over 45 random sam- 
ples). It is clear that both models lead to the same emission 
maps. While of course the actual intensity of the emitted TH 
signal in incoherent and coherent models will differ, in partic- 
ular for small dispersion case, the angular dependence of the 
emitted TH will be the same. Again we clearly see that the 
out of many discrete TH emission channels clearly visible 
in ideal periodic structure only the strongest, namely those 
involving Cerenkov and forward emission are preserved in 
highly disordered regime. 



5 Conclusions 

In conclusions, we investigated theoretically the second and 
third harmonic generation in periodically poled nonlinear crys- 
tals with random ferroelectric domain distribution. In partic- 
ular, we discuss properties of second and cascaded third har- 
monic generation. We derived analytical formulas describing 
emission properties of the second and third harmonics in the 
presence of domain disorder. We showed that in the limit of 
large number of domains, the THG process is described by 
product of simple expressions describing each of the con- 
stituent processes, i.e. SHG and sum frequency mixing. We 
considered various processes involved in the frequency mix- 
ing and analyzed their properties as a function of wavelength 
and degree of domain disorder. We show that the randomness- 
induced incoherence in the wave interaction leads to deteri- 
oration and angular spreading of harmonic generated in the 
processes relying on transverse phase matching such as Raman- 
Nath. On the other hand forward and Cerenkov frequency 
generation are basically insensitive to the domain random- 
ness. 
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